/*
 heat-spectral lib. solves time-dependent PDEs using SEM-NI method.
 Copyright (C) 2009, Altieri Michele and Stefanoni Andrea
 
 This program is free software: you can redistribute it and/or modify
 it under the terms of the GNU General Public License as published by
 the Free Software Foundation, either version 3 of the License, or
 (at your option) any later version.
 
 This program is distributed in the hope that it will be useful,
 but WITHOUT ANY WARRANTY; without even the implied warranty of
 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
 GNU General Public License for more details.
 
 You should have received a copy of the GNU General Public License
 along with this program. If not, see <http://www.gnu.org/licenses/>.
*/

#include<ctime>
#include<iostream>
#include<list>

#include "advdiff.hpp"

using namespace spectral;
using namespace std;

int main() {
    // =============================================
    // Problem Parameters
    // =============================================
    // Number of elements along the x-axis
    int nx(100);
    // Number of elements along the y-axis
    int ny(100);
    // Number of nodes
    int nnodes(4);
    // Time step
    double deltat(1./25);
    // Number of time steps
    int ntime(500);
    // Relaxation parameter
    double theta(1);
    // nmax for the newton method
    int newton_nmax((int)1e4);
    // Order of the problem in the week form
    int order(1);
    // Tolerance parameters for the newton method
    double legDelta(1e-15);
    double legEps(1e-15);
    // =============================================

    // =============================================
    // Solving the problem
    // =============================================
    clock_t ck0;
    clock_t ck1;

    cout << "Initializing SEM problem..." << flush;
    ck0=clock();
    Problem<ADV_PROBLEM> pb(LA,LB,nx,ny,nnodes,deltat,theta,newton_nmax,
			     legDelta,legEps,order);
    pb.initialize();
    pb.printSol("./advData/advSol0000.dat");

    ck1=clock();
    cout << "done! ";
    cout << "(In " << double(ck1-ck0)/CLOCKS_PER_SEC << " seconds)\n";
    cout << flush;

    cout << "Solving parabolic problem..." << endl;
    
    char filename[30];
    double actualtime(0);
    double solMin = pb.range(0);
    double solMax = pb.range(1);
    
    ck0=clock();
    
    for (int t=0;t<ntime;++t) {
	actualtime = deltat*(t+1);
	cout << "Solving the time " << actualtime << endl;
	pb.solve(actualtime);

	if (pb.range(0) < solMin) solMin = pb.range(0);
	if (pb.range(1) > solMax) solMax = pb.range(1);

	sprintf(filename, "%s%04d%s", "./advData/advSol", t+1, ".dat");
	pb.printSol(filename);
    }
    
    ck1=clock();
    cout << "done! ";
    cout << "(In " << double(ck1-ck0)/CLOCKS_PER_SEC << " seconds)\n";
    cout << flush;
    // =============================================

    // =============================================
    // Prints the scripts for Gnuplot and Mencoder
    // =============================================
    pb.printAviScript("adv","Diffusion-Advection Problem",ntime,solMin,solMax);
    pb.printGnuplotScript("adv","Diffustion-Advection Problem",solMin,solMax);
    // =============================================

    // =============================================
    // Results output
    // =============================================
    cout << endl;
    cout << "Minimum value of the solution: " << solMin << endl;
    cout << "Maximum value of the solution: " << solMax << endl;	       
    // =============================================
}
